{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Explore LAT Data (for Burst)\n",
    "\n",
    "In this example, we will examine the LAT data for a gamma-ray burst based on the time and position derived from a GBM trigger."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Prerequisites\n",
    "\n",
    "It is assumed that:\n",
    "\n",
    "* You are in your working directory.\n",
    "* The GBM reported a burst via a [GCN circular]() with the following information:\n",
    "    * Name = GRB 080916C\n",
    "    * RA = 121.8\n",
    "    * Dec = -61.3\n",
    "    * TStart = 243216766 s (Mission Elapsed Time)\n",
    "    * T90 = 66 s"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Steps\n",
    "\n",
    "The analysis steps are:\n",
    "\n",
    "1. Extract the Data\n",
    "2. Data Selections\n",
    "3. Bin the Data\n",
    "4. Look at the Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Extract the Data\n",
    "\n",
    "Generally, one should refer to [Extract LAT data](https://fermi.gsfc.nasa.gov/ssc/data/p6v11/analysis/scitools/extract_latdata.html) tutorial and use the time and spatial information from some source (here, a GCN notice) to make the appropiate extraction cuts from the [LAT data server](http://fermi.gsfc.nasa.gov/cgi-bin/ssc/LAT/LATDataQuery.cgi)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use the following parameters to extract the data for an ROI of 20 degrees from 500s before to 1500 seconds after the trigger time:\n",
    "\n",
    "* Search Center (RA,Dec)\t=\t(121.8,-61.3)\n",
    "* Radius\t=\t20 degrees\n",
    "* Start Time (MET)\t=\t243216266 seconds (2008-09-16T00:04:26)\n",
    "* Stop Time (MET)\t=\t243218266 seconds (2008-09-16T00:37:46)\n",
    "* Minimum Energy\t=\t20 MeV\n",
    "* Maximum Energy\t=\t300000 MeV\n",
    "\n",
    "Note that for analyses that require the diffuse background, the standard model starts at 60 MeV."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!wget https://fermi.gsfc.nasa.gov/FTP/fermi/data/lat/queries/L1906251643533A87413932_PH00.fits\n",
    "!wget https://fermi.gsfc.nasa.gov/FTP/fermi/data/lat/queries/L1906251643533A87413932_SC00.fits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!mkdir data\n",
    "!mv *.fits ./data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!mv ./data/*_PH00.fits ./data/grb_events.fits\n",
    "!mv ./data/*.fits ./data/grb_spacecraft.fits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!ls ./data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This will extract the data files used in this tutorial into the `data` directory."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data Selections\n",
    "\n",
    "NOTE: For information on the recommended selections for burst analysis of LAT data, you should refer to the [Cicerone](http://fermi.gsfc.nasa.gov/ssc/data/p6v11/analysis/documentation/Cicerone/Cicerone_Data_Exploration/Data_preparation.html).\n",
    "\n",
    "To map the region of the burst, you should select a large spatial region from within a short time range bracketing the burst. For the lightcurve on the other hand, its best to select you want a small region around the burst from within a long time range. Therefore, you will need to make two different data selections.\n",
    "\n",
    "In both cases, we will use the loosest event class cut (to include all class events). This works because rapid, bright events like GRBs overwhelm the background rates for the short period of the flare."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, we use [gtselect](https://fermi.gsfc.nasa.gov/ssc/data/p6v11/analysis/scitools/help/gtselect.txt) to extract the region for the spatial mapping:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "gtselect evclass=128\n",
    "    ./data/grb_events.fits\n",
    "    ./data/GRB081916C_map_events.fits\n",
    "    121.8\n",
    "    -61.3\n",
    "    20\n",
    "    243216666\n",
    "    243216966\n",
    "    100\n",
    "    300000\n",
    "    180"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Notes**:\n",
    "* The input photon file is `grb_events.fits`, and the output file is `GRB081916C_map_events.fits`.\n",
    "* We selected a circular region with radius 20 degree around the burst location (the ROI here has to fall within that selected in the data server) from a 300s time period around the trigger.\n",
    "* We made an additional energy cut, selecting only photons between 100 MeV and 300 GeV. This removes the low-energy high-background events from the map."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Run [gtselect](https://fermi.gsfc.nasa.gov/ssc/data/p6v11/analysis/scitools/help/gtselect.txt) again**, this time with a smaller region but a longer time range. (Note that _gtselect_ saves the values from the previous run and uses them as defaults for the next run.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "gtselect evclass=128\n",
    "    ./data/grb_events.fits\n",
    "    ./data/GRB081916C_lc_events.fits\n",
    "    121.8\n",
    "    -61.3\n",
    "    10\n",
    "    243216266\n",
    "    243218266\n",
    "    30\n",
    "    300000\n",
    "    180"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Notes**:\n",
    "\n",
    "* A new output file, with the extension `lc_events.fits` has been produced.\n",
    "* The search radius was reduced to 10 degrees.\n",
    "* The start to stop time range has been expanded.\n",
    "* The energy range has been expanded."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bin the Data\n",
    "\n",
    "Use [gtbin](https://fermi.gsfc.nasa.gov/ssc/data/p6v11/analysis/scitools/help/gtbin.txt) to bin the photon data into a map and a lightcurve.\n",
    "\n",
    "First, create the counts map:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "gtbin\n",
    "    CMAP\n",
    "    ./data/GRB081916C_map_events.fits\n",
    "    ./data/GRB081916C_counts_map.fits\n",
    "    NONE\n",
    "    50\n",
    "    50\n",
    "    0.5\n",
    "    CEL\n",
    "    121.8\n",
    "    -61.3\n",
    "    0.\n",
    "    AIT"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Notes**:\n",
    "    \n",
    "* Select _cmap_ (i.e., count map) as the output.\n",
    "* When we ran _gtselect_ we called the file from a large area `GRB081916C_map_events.fits`.\n",
    "* The counts map file will be called `GRB081916C_counts_map.fits`.\n",
    "* Although in _gtselect_ we selected a circular region with a 20 degree radius, we will form a counts map with 50 pixels on a side, and 1/2 degree square pixels.\n",
    "\n",
    "Now, we create the lightcurve. Once again, the previous inputs were saved and provided as defaults where appropriate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "gtbin\n",
    "    LC\n",
    "    ./data/GRB081916C_lc_events.fits\n",
    "    ./data/GRB081916C_light_curve.fits\n",
    "    NONE\n",
    "    LIN\n",
    "    243216266\n",
    "    243218266\n",
    "    10"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Notes**:\n",
    "* We used **gtselect** to create `GRB081916C_lc_events.fits`, a photon list from a small region but using a long time range.\n",
    "\n",
    "\n",
    "* **gtbin** binned these photons and output the result into `GRB081916C_light_curve.fits`.\n",
    "\n",
    "  There are a number of options for choosing the time bins. Here we have chosen linear bins with equal time widths of 10 seconds."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Examine the Data\n",
    "\n",
    "We now have two FITS files with binned data: one with a lightcurve (`GRB081916C_light_curve.fits`), and a second with a counts map (`GRB081916C_counts_map.fits`).\n",
    "\n",
    "**Note**: Currently we do not have any Fermi specific graphics programs, but there are various tools available to plot FITS data files, such as [_fv_](http://heasarc.nasa.gov/ftools/fv/) and [_ds9_](http://hea-www.harvard.edu/RD/ds9/)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To look at the counts map, we will use *fv*.\n",
    "\n",
    "1. First, start up *fv* in your terminal using:\n",
    "       prompt> fv &\n",
    "    (Or see the code cell below)\n",
    "       \n",
    "       \n",
    "2. Then open `GRB081916C_counts_map.fits`. Click on 'open file' to get the `File Dialog` GUI. Choose `GRB081916C_counts_map.fits`. A new GUI will open up with a table with two rows.\n",
    "\n",
    "\n",
    "3. FITS files consist of a series of extensions with data. Since the counts map is an image, it is stored in the primary extension (for historical reasons only images can be stored in the primary extension). Clicking on the 'Image' button in the first row results in:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!fv ./data/GRB081916C_counts_map.fits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from IPython.display import HTML"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "HTML(\"<img src='https://fermi.gsfc.nasa.gov/ssc/data/analysis/scitools/images/explore_data_GRB/GRB081916C_counts_map_fv.png'>\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One readily observes that there are many pixels containing a few counts each near the burst location, and few pixels containing any counts away from the burst. It is also evident that the burst location is not centered in the counts map. This means the preliminary location sent out in the GCN was not quite the right location. This is not uncommon for automated transient localizations.\n",
    "\n",
    "Using *ds9* instead of *fv* allows us to find the location of the brightest pixel, and a better position for the burst."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "HTML(\"<img src='https://fermi.gsfc.nasa.gov/ssc/data/analysis/scitools/images/explore_data_GRB/GRB081916C_counts_map_ds9.png'>\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By simply mousing over that pixel, the RA and Dec of the likely burst location get displayed in the FK5 boxes (here, the new position is 119.5, -56.5). This is a quick method of localizing a burst, though not as accurate as fitting the data (described in the [Likelihood Tutorial](https://fermi.gsfc.nasa.gov/ssc/data/p6v11/analysis/scitools/likelihood_tutorial.html)).\n",
    "\n",
    "> **Note**: To change the *ds9* display for the mouseover from sexagesimal (the default) to degrees, go to the WCS menu and select \"Degrees\" in the bottom portion of the menu."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**To look at the lightcurve:**\n",
    "\n",
    "1. Open `GRB081916C_light_curve.fits` in *fv*.\n",
    "\n",
    "    The lightcurve is in the RATE extension. Choose 'All' to view the content of that extension. The extension has 4 columns: TIME, TIMEDEL, COUNTS and ERROR.\n",
    "    \n",
    "    The time bins were created to have a width of 10 seconds, and therefore TIMEDEL is always 10.\n",
    "    \n",
    "2. Now plot the COUNTS as a function of TIME. Select 'Plot' for the RATE extension. Click on 'Time' then on 'x'. Click on 'Counts' then on 'y'. Finally click on 'Plot'. The result is:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!fv ./data/GrB081916C_light_curve.fits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "HTML(\"<img src='https://fermi.gsfc.nasa.gov/ssc/data/analysis/scitools/images/explore_data_GRB/GRB081916C_light_curve_fv.png'>\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    ">**Note**: Except for the period of the burst, almost all bins have 0-5 counts. This demonstrates that there is very little background for LAT observations of gamma-ray bursts."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
